function pyr = pyramid_sim(base, filt_coeff, levels, gaus_or_lap)
% gaus_or_lap  -  0 : gaussian (default), 1 : laplacian

dec_factor = 2; % dyadic pyramid assumed

if(gaus_or_lap == 0)    %for gaussian pyramid

    pyr{1} = base;

    for lev = 2:levels+1

        pyr_lev = pyr{lev-1};
        
        filter_out = conv2_sim(pyr_lev, filt_coeff); %CONV2
                
%         filter_out = imfilter_sim(pyr_lev, filt_coeff);   %IMFILTER 
        
        lpdxy_out = decimate_sim(filter_out, dec_factor, 0);    
        %pixel_select - 0 : select first pixel in the window and drop the
        %rest, 1 : select last pixel in the window and drop the rest
        pyr{lev} = lpdxy_out;

    end

else                       % 1 for laplacian pyramid
    
    
    gpyr{1} = base;

    for lev = 2:levels+2

        pyr_lev = gpyr{lev-1};
        filter_out = conv2_sim(pyr_lev, filt_coeff);    %CONV2
        
%         filter_out = imfilter_sim(pyr_lev, filt_coeff);   %IMFILTER 
        
        pyr {lev-1} = pyr_lev - filter_out;    %lap(i) = gauss(i) - conv2(gauss(i))
        
        lpdxy_out = decimate_sim(filter_out, dec_factor, 0);
        %pixel_select - 0 : select first pixel in the window and drop the
        %rest, 1 : select last pixel in the window and drop the rest
        gpyr{lev} = lpdxy_out;

    end
  
end
    